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ABSTRACT 

A three-dimensional flow analysis code has been 
used to compute the design speed operating line of a 
transonic fan rotor, and the results have been com- 
pared with experimental data. The code is an explicit 
finite difference code with an algebraic turbulence 
model. It is described in detail in the paper. The tran- 
sonic fan, designated rotor 67, has been tested exper- 
imentally at NASA Lewis Research Center with con- 
ventional aerodynamic probes and with laser anemom- 
etry and has been included as one of the AGARD 
test cases for the computation of internal flows. The 
experimental data are described briefly. Maps of to- 
tal pressure ratio and adiabatic efficiency versus mass 
flow have been computed and are compared with the 
experimental maps, with excelent agreement between 
the two. Detailed comparisons between calculations 
and experiment are made at two operating points, one 
near peak efficiency and the other near stall. Blade-to- 
blade contour plots are used to show the shock struc- 
ture. Comparisons of circumferentially integrated flow 
quantities downstream of the rotor show spanwise dis- 
tributions of several aerodynamic parameters. Calcu- 
lated Mach number distributions are compared with 
laser anemometer data within the blade row and the 
wake to quantify the accuracy of the calculations. Fi- 
nally, particle traces are used to illustrate the nature 
of the secondary flow in this fan. 

INTRODUCTION 

One goal of computational fluid dynamics (CFD) 
for turbomachinery is the prediction of component per- 
formance, for example pressure ratio and efficiency. 
Since a small improvement in engine efficiency can 
amount to huge savings in yearly fuel costs for a fleet of 
commercial aircraft, turbomachinery designers are ex- 
tremely interested in tools that give good quantitative 
predictions of turbomachinery performance. However, 
most of the CFD results for turbomachinery published 
more than a few years ago show only qualitative com- 
parisons with experimental Mach contours or quantity 
tive comparisons with surface pressures. This is not to 
diminish the importance of these results; in the hands 
of a good turbomachinery designer a shock location or 
surface pressure distribution can be extremely useful. 
However, CFD for turbomachinery has evolved to a 


point where it can now provide reasonable predictions 
of overall performance of isolated blade rows. 

In more recent literature several researchers have 
presented more detailed predictions of turbomachin- 
ery performance. Davis et al. (ref. 1) predicted loss 
buckets for two-dimensional transonic compressor cas- 
cades. Chevrin and Vuillez (ref. 2) have predicted loss 
and exit flow angle for turbine and fan cascades using 
the two-dimensional code of Cambier et al. (ref. 3). 
Boyle (ref. 4) investigated effects of turbulence mod- 
eling on turbine blade heat transfer predictions using 
the two-dimensional code of Chima (ref. 5). In three 
dimensions, Pierzga and Wood (ref. 6) predicted the 
static pressure ratio versus mass flow curve of a tran- 
sonic fan using Denton’s code (ref. 7), and Dawes (ref. 
8) has predicted exit total pressure and temperature 
distributions in a multistage turbine. Adamczyk, et 
al. (ref. 9) investigated the effects of tip clearance on 
stall for the fan considered in the present work and 
predicted pressure ratios and efficiencies for that fan. 
In reference 10 Chima used the present code to predict 
the overall efficiency of an annular turbine stator. 

There are many reasons for the scarcity of turbo- 
machinery performance calculations in the early liter- 
ature. One is that pressure field calculations are rela- 
tively independent of viscous effects and can be calcu- 
lated with simple models. Conversely, efficiency and 
loss calculations are highly dependent on viscous ef- 
fects and require careful attention to the viscous terms, 
turbulence modeling, artificial viscosity, and grid res- 
olution for successful calculation. Second, computers 
capable of performing large viscous flow calculations 
have become generally available only recently. A final 
reason is that detailed experimental data are difficult 
to obtain in turbomachinery because of the small size 
and high speeds of the components involved. Experi- 
mental surface pressures are available from many lin- 
ear and annular cascade tests, which may account for 
many of the comparisons appearing in the literature. 
However, wake surveys and loss data are often avail- 
able for the same tests, but tend to be overlooked. 

With the publication of AGARD Advisory Report 
No. 275, Test Cases for Computation of Internal Flows 
in Aero Engine Components (ref. 11), researchers now 
have access to some excellent data for the validation 
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of CFD codes for turbomachinery. 

In the present work the experimental data de- 
scribed by Strazisar et al. (ref. 12) and by Wood 
et al. (ref. 11, pp. 165-213) was used to validate the 
three-dimensional Navier-Stokes code first referenced 
by Chima and Yokota in (ref. 10). The code, EVC3D 
(Rotor Viscous Code 3D), is described with emphasis 
on the boundary conditions and artificial viscosity. A 
new three-dimensional grid code for turbomachinery 
is also introduced. The test case, a transonic fan rotor 
(rotor 67, shown in fig. 1) is described briefly. 

Several operating points were computed along the 
100 percent speed line of the rotor. The computed op- 
erating curves of adiabatic efficiency and total pressure 
ratio versus mass flow are compared with the exper- 
imental data. Two operating points are examined in 
detail: one near peak efficiency and one near stall. At 
each point qualitative comparisons are made between 
computed and experimental Mach number contours. 
Comparisons are then made with laser anemometer 
measurements within the blade row and across the 
wake, and with conventional aerodynamic measure- 
ments downstream. Finally, particle traces are used 
to illustrate the nature of the secondary flows in this 
machine. 

GOVERNING EQUATIONS 

The Navier-Stokes equations are written in a 
Cartesian (x, y, z) coordinate system rotating with an- 
gular velocity Q, about the x-axis. The rotation intro- 
duces source terms in the y and z momentum equa- 
tions. The Cartesian equations are mapped to a gen- 
eral body-fitted (£,*7,0 coordinate system using stan- 
dard techniques. A C-type grid was used in the present 
work, with the ^-coordinate roughly following the flow, 
the ^coordinate running blade-to-blade, and the £- 
coordinate running spanwise. The thin-layer approxi- 
mation is used to drop all viscous derivatives in the £- 
direction, All viscous terms in the cross-channel (f?,C) 
plane are retained. The resulting equations are as fol- 
lows: 

d t q + J[d<:E + d n F + d < G-Re- 1 (d r ,Fv + d< : Gv)] = H 

( 1 ) 
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Again, these equations are written in a rotating 
Cartesian coordinate system attached to the moving 
blade. The code solves for the absolute velocity com- 
ponents u,v, and w y which point in the x,y, and z 
coordinate directions, respectively. The relative veloc- 
ity components u' , t/, and w f are defined with respect 
to the same rotating coordinate system by subtracting 
the appropriate components of the blade speed from 
the absolute velocities, giving 


v ' — v — Qz (3) 

w i = w + Qy 


The relative contravariant velocity components 
are given by 


U* - + £yt/ +6™' 

V' .3S T) x u -h T} y v' + J] z uf (4) 

W' + Cy"' +<zW* 


Note that although u* — u y U* ^ U . The energy and 
static pressure are given by 
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Using Stokes’ hypothesis, A = — §p, the viscous 
flux Fy can be written as follows: 


F v = J'VfO, F 2 , F 3 , F 4 , F 5 ] t (7) 
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Terms multiplied by C\ and C 2 lead to nonmixed, sec- 
ond derivative viscous terms like u nn , while terms mul- 
tiplied by C 3 -C 5 lead to mixed-derivative terms like 
u v {. The viscous flux vector G can be written simi- 
larly by interchanging directions r/ and C and replacing 
F with G everywhere. 

Metric terms are defined using the following rela- 
tions: 
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Terms like etc. are found using second order 

central differences. The Jacobian is computed using 
equation (11) and stored for the entire grid. All other 
metric terms are computed as needed using equation 

(io). 

The equations are nondimensionalized by arbi- 
trary reference quantities (here, the inlet total density 
Poref ^d the total sonic velocity c 0re / were used), and 
the Reynolds number Re and Prandtl number Pr are 
defined in terms of these quantities. The equations as- 
sume that the specific heats C p and C v and Prandtl 
number are constant, that Stoked hypothesis is valid, 
and that the effective viscosity for turbulent flows may 
be written as 


Mi = Mm + Hturb ( 12 ) 

where the laminar viscosity is calculated using a power 
law function of temperature: 

(i3) 

with n = 2/3 for air. The turbulent viscosity pturb is 
computed using an adaptation of the Baldwin-Lomax 
turbulence model (ref. 13) on cross-channel planes. 
Briefly, the model is applied independently in the 
blade-to-blade and spanwise directions, and the re- 
sulting turbulent viscosities are added vectorally. The 
two-dimensional Buleev length scale based on the dis- 
tance from the hub and blade is used. The vorticity 
is calculated in the absolute frame, but the wall shear 
and wake velocities are calculated in the relative frame. 
See reference 10 for details on the three-dimensional 
implementation. 

BOUNDARY CONDITIONS 

Many boundary conditions for turbomachinery 
are best expressed in cylindrical coordinates, but the 


code is formulated in Cartesian coordinates. For the 
boundary conditions, coordinate systems are trans- 
formed as needed using: 

v = ( v r y+v $ z)/r 

w = (v r z-v 6 y)/r 
v r - ( vy + wz)/r (14) 

v e — ( vz - wy)/r 

r = \Jy 2 + z 2 

At the inflow boundary the total temperature, 
total pressure, whirl, and meridional flow angle were 
specified, and the upstream-running Riemann invari- 
ant was extrapolated from the interior. The inlet total 
temperature To was specified as a constant (standard 
conditions). The inlet total pressure was specified as 
a constant in the core flow (also standard conditions) 
and reduced in the endwall regions according to a 1/7 
power law velocity profile, with the inlet boundary 
layer heights estimated from the experimental data 
to be 12 mm on both the hub and the casing. At 
the inlet v$ was set to zero, and v r was chosen to 
make the flow tangent to the meridional projection 
of the inlet grid lines. The upstream-running Rie- 
mann invari ant R~ based on the total absolute ve- 
locity Q — y/u 2 -I- v 2 + w 2 was calculated at the first 
interior point and extrapolated to the inlet. The Rie- 
mann invariant is given by 

R-=Q-JL. (15) 

The total velocity is found from To and R~ using 

n (T — 1)-R~ + \/2(l — 7)(i?~) 2 4- 4(7 + 1 )C p Tq 
Q 7+1 

(16) 

The velocity components are then decoupled alge- 
braically, and the density is found from Po ) T 0) and 
Q using an isentropic relation. 

At the exit the hub static pressure is specified, and 
p ) pu y pv y and pw are extrapolated from the interior. 
The local static pressure is found by integrating the 
axisymmetric radial equilibrium equation: 


dp 

dr 


pve _ P f \ 2 

= -~{vz - wy) 


(17) 


A periodic C-grid was used in the present work. 
The periodic boundary was solved by setting periodic 
flow conditions (in terms of cylindrical velocity com- 
ponents) on a dummy grid line outside the boundary. 
Fourth-difference artificial dissipation terms are ne- 
glected on the outer boundary so that only one dummy 
grid line is needed. 

On the blade surface and the rotating part of the 
hub, U f — V f — W f — 0. The hub was specified to be 
stationary 13.4 mm upstream of the leading edge and 
3.35 mm downstream of the trailing edge by setting 
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u = v = w = 0. The tip casing was assumed to be sta- 
tionary, that is, the blade was assumed to scrape along 
the casing with zero tip clearance. Surface pressures 
were found from the normal momentum equation. On 
the hub « = 1) and tip (( = (mar) 

(Cr(r + (y(y T C t£z)dtP + (Cr^x + (y% + 

+ (Cr + Cy + Cz)d(P = -p[Q(( y W - Cv)] (18) 

On the blades (77 = 1) the normal momentum 
equation can be found by replacing ( with 77 every- 
where in equation (18). 

Surface temperatures were found from an adia- 
batic wall condition implemented as in equation (18) 
with p replaced by T and a right-hand side of zero. 

MULTISTAGE RUNGE-KUTTA SCHEME 

The governing equations are discretized using a 
node-centered finite difference scheme. Second-order 
central differences are used throughout. The mul- 
tistage Runge-Kutta scheme developed by Jameson, 
Schmidt, and Turkel (ref. 14) is used to advance the 
flow equations in time from an initial guess to a steady 
state. If equation (1) is rewritten as 

d t q = -J[Ri-(Rv + D)] (19) 

where Rj is the inviscid residual including the source 
term, Ry is the viscous residual, and D is an artificial 
dissipation term described in the next section, then 
the multistage Runge-Kutta algorithm can be written 
as follows: 

<IQ — ( ln 

qi = <70 — J At [Rj qo — ( Rv + D) go] 

: (20) 

qk = 90 - Rk 

= qo — akJ At [Ri gjb-i — ( Rv + D) qo] 
qn + 1 = qk 

Here a standard four-stage scheme was used, with 
a,- = 1/4, 1/3, 1/2, 1. For efficiency both the physical 
and artificial dissipation terms are calculated only at 
the first stage, then are held constant for subsequent 
stages. 

ARTIFICIAL DISSIPATION 

The dissipative term D in equation (19) is similar 
to that used by Jameson et al. (ref. 14). It is given by 

Dq = (D ( +D tl +D ( )q (21) 

where the ^-direction operator is given by 

D ( q = Q {V 3 q« - V 4 q^) (22) 

The terms Vi and V 4 are given by 

V 2 = PI +/i 2 max(t/ i+ i,i'i,J/i_i) 

V 4 — max (0, p 4 — V2) (23) 


where 

Ui = I^. + I -2P. + P-!! 

mm(P in ,\P i+ i+2P i + P i _i\) V ’ 

and subscript i corresponds to the ( direction, etc. 

The constant p,\ scales a first-order artificial 
viscosity that is useful for maintaining stability at 
startup. The constant m scales a first-order artifi- 
cial artificial viscosity that is switched on at shocks 
detected by equation (24). The denominator in equa- 
tion (24) is normally constant at the inlet pressure 
Pj n , making the operator roughly symmetric across 
shocks. The more common term |P,-fi T 2 P, -j- Pj_i| 
is included to switch on the second-difference dissipa- 
tion when the pressure becomes very small, usually 
due to numerical problems. The constant p 4 scales a 
uniform third-order artificial viscosity that is switched 
off at shocks by equation (23). In this work was set 
to 1/4 for the first 200 iterations and zero thereafter, 
pi = 1/2, and = 1/32. 

C is an arbitrary coefficient that can have a large 
impact on the stability and accuracy of the solution. 
The subscript indicates that C may be different in each 
direction. 

In reference 10 a directionally homogeneous coef- 
ficient was used for three-dimensional problems. To 
minimize dissipation in viscous regions C was reduced 
to zero linearly over sevaral grid points near walls. In 
reference 5 a directionally biased coefficient was used 
for two-dimensional problems. The dissipation was 
proportional to the grid spacing in each direction, thus 
reducing it across finely gridded viscous regions. This 
worked well in two dimensions but did not generalize 
to three dimensions. 

Martinelli and Jameson (ref. 15) proposed a direc- 
tionally biased coefficient that works well but can still 
have large directional variations on highly stretched 
three-dimensional grids. Kunz and Lakshminarayana 
(ref. 16) proposed a modified form of the Martinelli 
and Jameson coefficient: 



etc., where cr = 2/3. 

In the present work the one- dimensional time 
steps are written as follows: 

A*e= (veySTiFfiF)’ 1 (26) 

etc, where the normalized velocity scale is approx- 
imated as one. Dissipation coefficients given by equar 
tions (25) and (26) seem to give a good distribution of 
the dissipative terms in each direction and have been 
used throughout this work. 
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THREE-DIMENSIONAL STABILITY LIMIT 

The following expression is used for the time step: 

A* 

At < 

/*M + iy\v\ + />'l + y/c*(i* + i* + i*) + n* 

where 

lx = |£x| + |*?x| + iCri 

ly = l£y I i + iCy I (27) 

1 2 = If* I + |^| + IGI 

and A* is the maximum Courant number for the par- 
ticular multistage scheme. For the standard four-stage 
scheme A* « 2.8. 

To accelerate convergence to a steady state, the 
maximum permissible time step at each grid point was 
used giving a constant Courant number everywhere. 
The time step was updated every 50 iterations. 

IMPLICIT RESIDUAL SMOOTHING 

To further accelerate convergence it is desirable 
to use a time step even larger than the stability limit 
given by equation (27). To maintain stability, the 
residual calculated in equation (19) is smoothed af- 
ter each Runge-Kutta stage by an implicit smoothing 
operator, that is, 

(1 — e$6^)( 1 — e f? 6 r?r? )(l — e^8^)Rk = Rk (28) 

where 6 <5^, and <5^ are standard second difference 
operators and and e ^ are smoothing parameters. 

Linear stability analysis shows that the Runge- 
Kutta scheme may be made unconditionally stable us- 
ing implicit residual smoothing if the smoothing pa- 
rameters € are made sufficiently large (ref. 17). In one 
dimension 



gives unconditional stability if A* is the Courant limit 
of the unsmoothed scheme, and A is a larger operating 
Courant number. In three dimensions different e’s may 
be used in each direction, and their magnitudes may 
often be reduced below the value given by equation 
(29). Courant numbers, A* « 2.8 and A = 5.0, and 
smoothing parameters, = 0.4, e n = 0.55, and — 

0.45, were used in the present work. 

COMPUTATIONAL GRID 

A three-dimensional grid code for turbomachinery 
has been developed by the author (currently unpub- 
lished.) The code, called TCGRID for Turbomachin- 
ery C-GRID, generates three-dimensional C-type grids 
using the following technique: 

1. A coarse, equally spaced meridional grid is 
generated between the specified hub and tip. 

2. Blade coordinates are found at the meridional 
grid points by interpolation of the input blade geome- 
try. 


3. Two-dimensional blade- to- blade grids are gen- 
erated along the meridional grid lines in ( m,r# ) co- 
ordinates using a version of the GRAPE code devel- 
oped by Steger and Sorenson (ref. 18). Here, m is 
the arc length along the meridional surface, and r is 
some mean radius. The GRAPE code allows arbitrary 
specification of inner and outer boundary points, then 
generates interior points as the solution of a Poisson 
equation. Forcing terms in the Poisson equation are 
chosen to maintain the desired grid spacing and angles 
at the boundaries. 

4. The (m, rd) coordinates are transformed back 
to (z,r,0). 

5. The two-dimensional grids are reclustered 
spanwise using a hyperbolic tangent stretching func- 
tion to make a full three-dimensional grid. 

6. Finally, the (z,r,0) coordinates are trans- 
formed to (x,y, z) and stored in a standard format. 

Figure 2 shows the 185x40x49 grid used for rotor 
67. The C-shaped grids used here give good resolution 
of the round leading edge of the blade, as shown in 
figure 3. The initial grid spacing is about 0.015 mm at 
the blade, 0.03 mm at the hub and 0.045 mm at the 
tip. 

COMPUTATIONAL DETAILS 

All computations were run on the Cray Y-MP 
computer at NASA Ames Research Center, under 
support from the Numerical Aerodynamic Simulation 
(NAS) Project Office. 

The grid code required about 3 million words 
(Mw) of in-core storage and ran in about 15 seconds 
for the grid shown in figures 2 and 3 (362 600 points.) 
The flow solver required about 6.5 Mw of storage and 
3.6 Mw of solid state device storage (SSD.) The SSD 
storage is used to hold the old solution go and the dis- 
sipative terms Ry + D during the four stages of the 
multistage scheme (eq. (20)). 

The spanwise inlet profile was used as an initial 
guess, and 250 iterations were run with an exit pres- 
sure corresponding to the peak efficiency point. About 
20 minutes of CPU time were required for this startup 
solution. All subsequent calculations were restarted 
from this solution and run an additional 1550 itera- 
tions, requiring about 2.5 CPU hours per case. 

EXPERIMENTAL DETAILS 

Experimental details are described briefly below. 
Full details may be found in references 11 and 12. 

Test rotor. The test rotor, NASA designation ro- 
tor 67, is shown in figure 1. It is the first stage rotor of 
a two stage fan, with a design pressure ratio of 1.63 at 
a mass flow of 33.25 kg/sec. The rotor has 22 blades. 
The tip radius varies from 25.7 cm at the leading edge 
to 24.25 cm at the trailing edge, and the hub to tip 
radius ratio varies from 0.375 to 0.478. At the design 
rotational speed of 16 043 rpm the tip speed is 429 
m/sec, and the tip relative Mach number is 1.38. 
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Aerodynamic performance measurements. The ro- 
tor mass flow was determined using a calibrated ori- 
fice. Radial surveys of total pressure and temperature, 
static pressure, and flow angle were made 2.54 cm up- 
stream of the leading edge and 2.02 cm downstream of 
the trailing edge using conventional probes. 

Laser anemomeiry measurements . A single chan- 
nel fringe anemometer was used for the measurements 
described below. Optical access to the fan was through 
a glass window in the casing. Fluorescent oil droplets 
with a diameter of about 1.0-1. 4 pm were used to seed 
the flow. 

Measurements were made by holding the probe 
volume fixed in space as the blades rotated by. Mea- 
surements were taken in 50 intervals going from the 
suction surface of one blade to the suction surface of 
the next blade. In references 11 and 12, and in the 
present work the term windows is used interchangeably 
with the term intervals in referring to laser measure- 
ment locations. Within the blade row the blade itself 
obscures the last few windows. Measurements were 
taken in 17 consecutive blade passages, then averaged. 
Only axial and tangential velocities were measured; 
however, Pierzga and Wood (ref. 6) have shown nu- 
merically that neglecting the radial component has a 
minimal effect on calculated Mach numbers. 

RESULTS 

The computed total pressure ratio and adiabatic 
efficiency at 100 percent speed are plotted against nor- 
malized mass flow and compared with experimental 
data in figure 4. The computed and experimental 
mass flows are normalized by their respective chok- 
ing mass flows, as suggested by Pierzga and Wood. 
This normalization removes any uncertainties in the 
experimental mass flows. The choking mass flow was 
measured as 34.96 kg/sec using a calibrated orifice, 
and the computed value was 34.54 kg/sec, a difference 
of 1.2 percent. The computed pressure ratios and ef- 
ficiencies agree very well with the experimental data, 
except that they are slightly high near stall. 

Figures 5 to 8 show detailed results at an operat- 
ing point near peak efficiency. Figure 5 shows a com- 
parison of relative Mach number contours at 10, 30, 
and 70 percent span from the tip. The experimental 
contours were drawn from laser anemometry data us- 
ing smoothing and interpolation procedures outlined 
in reference 6. At 10 percent span the inlet Mach 
number is about 1.35, and a bow wave stands ahead 
of the blade. A weak oblique shock crosses the pas- 
sage inside the blade row, and a strong normal shock 
sits near the trailing edge. The flow exits at a Mach 
number of about 0.95. The flow is qualitatively simi- 
lar at 30 percent span, except that the shock sits more 
forward. At 70 percent span the inlet Mach number 
is about 0.95. A small supersonic bubble forms on the 
forward portion of the suction surface. It is not clear 
if this bubble is terminated by a shock. 


Figure 6 shows a comparison of computed and 
measured blade-to-blade profiles of relative Mach num- 
ber versus laser window number at two axial locations 
for each of three spanwise locations. The first axial lo- 
cation is within the blade row and the second is within 
the wake. Again, 50 laser windows are defined from 
the suction surface of one blade (window 0) to the suc- 
tion surface of the next blade (window 50.) 

Within the blade row the computed profiles agree 
reasonably well with the laser data and predict the 
shock locations accurately. At 70 percent span near 
the pressure surface (laser windows greater than 40), 
few seed particles were observed, and the data are sta- 
tistically uncertain. 

Within the wake the computed Mach numbers are 
somewhat low near the tip and get progressively bet- 
ter towards the hub. The computed wake profiles are 
deeper than the measured profiles. It is thought that 
this is because the seed particles used for the laser 
anemometer measurements could not follow the high 
shear rates found at the center and edges of the wake. 

The computational results were averaged blade- 
to-blade using an “energy average” procedure devel- 
oped by D. L. Tweedt at NASA Lewis Reasearch Cen- 
ter (unpublished.) On each blade-to-blade grid line 
the procedure integrates the mass flow, radial and tan- 
gential momentum, total enthalpy, and ideal total en- 
thalpy P 0 y . The radial distributions were mass aver- 
aged spanwise to produce the overall averages shown 
previously in figure 4. Since the integrations give total 
conditions directly, they can be expected to give good 
predictions of efficiency. Other quantities like aver- 
age static pressure or flow angle are found as nonlin- 
ear algebraic combinations of the primitive integrated 
quantities, and thus may not agree as well with exper- 
imental data. 

Experimentally, the exit total conditions and flow 
angle were measured with a self-nulling combination 
probe, and static pressure was measured using a sep- 
arate self-nulling wedge probe. Total conditions are 
thought to be fairly accurate, but flow angle and static 
pressure accuracy depend on the frequency response 
and blockage of the instrumentation. 

Radial surveys of several aerodynamic parameters 
measured 2.02 cm downstream of the rotor (at the grid 
exit) are shown in figure 7. The computed exit total 
temperatures and total pressures agree very well with 
the measurements along the span, but the computed 
static pressures are somewhat high, and computed exit 
flow angles are two to three degrees high over much of 
the span. 

Figure 8 shows several particle traces to illustrate 
some of the secondary flows in this fan. The overall 
view shows particles introduced upstream so as to pass 
over the suction surface of the blade. Most of the par- 
ticles pass straight through the blade row, but those in 
the endwall boundary layer roll up into a vortex that 
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climbs the leading edge. The enlargement of the lead- 
ing edge shows that the flow in the hub boundary layer 
sees a high blade incidence and separates near the lead- 
ing edge. The low-momentum fluid in this separated 
region is centrifuged radially outwards. Eventually the 
incidence decreases, and the flow turns abruptly in the 
streamwise direction. The trailing edge enlargement 
shows a large separation bubble fed by fluid from the 
endwall that migrates radially outward and ends up 
in the wake. The complicated flow pattern in this re- 
gion may explain the scarcity of laser seed particles 
observed in this region. 

Figures 9 to 12 show detailed results at an oper- 
ating point near stall. Figure 9 shows a comparison 
of relative Mach number contours, again at 10, 30, 
and 70 percent span from the tip. At 10 percent span 
the inlet Mach number is about 1.4. A normal shock 
stands ahead of the blade and crosses the passage. 
The computed shock is somewhat stronger than, and 
ahead of, the measured shock. The exit Mach num- 
ber is about 0.85. At 30 percent span the results are 
similar. The flow at 70 percent span is similar to the 
peak efficiency case, except that the supersonic bubble 
is smaller and has a better defined (though smeared) 
terminating shock. 

Figure 10 shows a comparison of relative Mach 
number profiles, again at two axial locations and three 
spanwise locations. At the tip the computed shocks are 
ahead of and stronger than the measured shocks, which 
accounts for the high predicted pressure ratios near 
stall. Adamczyk et al. (ref. 9) have shown numerically 
that interaction of the tip leakage vortex and the tip 
shock has a large effect on the near-stall performance 
of this rotor, and that lack of a tip clearance model 
can account for the discrepancies in shock position and 
strength seen here. At 70 percent the computed Mach 
numbers are slightly low, but agree qualitatively with 
the data. 

In the wake, the computed Mach numbers are 
somewhat low near the tip and get progressively bet- 
ter towards the hub. Again, the computed wakes are 
deeper than the measured wakes. There is consid- 
erable uncertainty in the laser data in the center of 
the wake at 70 percent span, where computed particle 
traces shown later indicate a large separation. 

Figure 11 compares radial surveys downstream of 
the rotor. The computed total temperatures and pres- 
sures are slightly high along the span, consistent with 
the high adiabatic efficiency shown in figure 4. Again, 
the computed static pressures are high. Computed 
exit flow angles agree well near the tip but are a few 
degrees high at lower radii. 

Figure 12 shows particle traces at the near-stall 
operating point. The overall view shows more radial 
migration of the endwall flow than at peak efficiency. 
At this lower mass flow the blade sees a higher relative 
incidence than at peak efficiency, causing the flow near 
the hub to migrate tangentially away from the suc- 


tion surface, as shown in the leading edge enlargement. 
The trailing edge enlargement shows a separation bub- 
ble similar to that seen near peak efficiency. The trail- 
ing edge separation appears to be slightly larger at the 
near-stall point. 

SUMMARY 

A flow analysis code has been developed for three- 
dimensional viscous flows in turbomachinery. The 
analysis solves the Navier-Stokes equations written in 
a general body-fitted coordinate system, including ro- 
tation about the ar-axis. The thin-layer approxima- 
tion is made in the streamwise direction but all vis- 
cous terms are included in the cross planes. The 
Baldwin-Lomax eddy-viscosity model is used for tur- 
bulent flows. 

An explicit multistage Runge-Kutta scheme is 
used to solve the finite- difference form of the flow 
equations. A variable time step and implicit resid- 
ual smoothing are used to accelerate the convergence 
of the scheme. The code is highly vectorized for the 
Cray Y-MP, and solutions can be computed on fairly 
fine grids in 2 to 3 hours. 

The code was used to compute the operating map 
of a transonic fan at design speed and showed good 
agreement with measured values of total pressure ra- 
tio and adiabatic efficiency. The computed results are 
slightly optimistic near stall, probably due to the lack 
of a tip clearance model. 

Detailed comparisons were made with experimen- 
tal data at two operating points, one near peak effi- 
ciency and one near stall. Comparisons were made 
with aerodynamic surveys downstream of the fan. In 
general, exit total temperature and total pressure were 
predicted quite accurately, but static pressure and flow 
angle showed some disagreement with the data. Com- 
parisons were also made with laser anemometry data 
in the blade row and in the wake. Shock location and 
strength were predicted closely near peak efficiency, 
but were overpredicted near stall. Predicted wake pro- 
files had about the right location and spread but were 
much deeper than measured wakes, probably due to 
lack of resolution in the laser data. 

Particle traces showed separated flow at both the 
leading and trailing edges at both operating points. 
The leading edge separation is too small to be seen 
in the laser data. The trailing edge separation may 
have been suggested indirectly by the scarcity of seed 
particles in this region. The ability of the code to 
predict these separated flow features suggests that the 
code could be used to guide experimental work aimed 
at resolving these features or to eliminate such features 
during design. 

Overall the code showed very good agreement 
with a variety of experimental data, thereby increas- 
ing confidence that the code can reliably be used to 
predict the performance of other machines as well. 


7 



References 

[1] Davis, R. L., Hobbs, D. E., and Weingold, H. D., 
“Prediction of Compressor Cascade Preformance 
Using a Navier-Stokes Technique,” ASME Paper 
88-GT-96, June 1986. 

[2] Chevrin, P.-A., and Vuillez, C., “Viscous 
Flow Computations in Turbomachine Cascades,” 
ASME Paper 90-GT-76, June 1990. 

[3] Cambier, L., Escande, B., and Veuillot, J. P., 
“Computation of Internal Flows at High Reynolds 
Number by Numerical Solution of the Navier- 
Stokes Equations,” Reel i. Aerosp. (English ed.), 
No. 1986, pp. 27-44. 

[4] Boyle, R. J., “Navier-Stokes Analysis of Turbine 
Blade Heat Transfer,” ASME Paper 90-GT-42, 
June 1990. 

[5] Chima, R. V., “Explicit Multigrid Algorithm for 
Quasi-Three-Dimensional Viscous Flows in Tur- 
bomachinery,” 7. Propulsion Power , Vol. 3, No. 
5, Sept.-Oct. 1987, pp. 397-405. 

[6] Pierzga, M. J., and Wood, J. R., “Investigation of 
the Three-Dimensional Flow Field Within a Tran- 
sonic Fan Rotor: Experiment and Analysis,” 7. 
Eng . Gas Turbines Power , Vol. 107, No. 2, Apr. 
1985, pp. 436-449. 

[7] Denton, J. D., “An Improved Time Marching 
Method For Turbomachinery Flow Calculations,” 
ASME Paper 82-GT-239, Apr. 1982. 

[8] Dawes, W. N., “Towards Improved Throughflow 
Capability: The Use of 3-D Viscous Flow Solvers 
in a Multistage Environment,” ASME Paper 90- 
GT-18, June 1990. 

[9] Adamczyk, J. J., Celestina, M. L., and Greitzer, 
E. M.*, “The Role of Tip Clearance in High-Speed 
Fan Stall,” ASME Paper 91-GT-83, June 1991. 

[10] Chima, R. V., and Yokota, J. W., “Numerical 
Analysis of Three-Dimensional Viscous Internal 
Flows,” AIAA 7., Vol. 28, No. 5, May 1990, pp. 
798-806. 

[11] Fottner, L., ed., Test Cases for Computation 
of Internal Flows in Aero Engine Components , 
AGARD Advisory Report No. 275, AGARD, 
Neuilly-Sur-Seine, France, July 1990. 

[12] Strazisar, A. J., Wood, J. R., Hathaway, M. D., 
Suder, K. L., “Laser Anemometer Measurements 
in a Transonic Axial Flow Fan Rotor,” NASA 
TP-2879, Nov. 1989. 

[13] Baldwin, B. S., and Lomax, H., “Thin-Layer 
Approximation and Algebraic Model for Sepa- 
rated Turbulent Flows,” AIAA Paper 78-257, 
Jan. 1978. 

[14] Jameson, A., Schmidt, W., and Turkel, E., “Nu- 


merical Solutions of the Euler Equations by Fi- 
nite Volume Methods Using Runge-Kutta Time- 
Stepping Schemes,” AIAA Paper 81-1259, June 
1981. 

[15] Martinelli, L., and Jameson, A., “Validation of 
a Multigrid Method for the Reynolds Averaged 
Equations,” AIAA Paper 88-0414, Jan. 1988. 

[16] Kunz, R. F., and Lakshminarayana, B. “Compu- 
tation of Supersonic and Low Subsonic Cascade 
Flows Using an Explicit Navier-Stokes Technique 
and the k — e Turbulence Model,” in Computa- 
tional Fluid Dynamics Symposium on Aeropropul- 
sion , NASA CP-10045, Apr. 1990. 

[17] Jameson, A., and Baker, T. J., “Solution of the 
Euler Equations for Complex Configurations,” 
AIAA Paper 83-1929, July 1983. 

[18] Steger, J. L., and Sorenson, R. L., “Auto- 
matic Mesh-Point Clustering Near a Boundary in 
Grid Generation with Elliptic Partial Differential 
Equations,” 7. Comp. Phys. f Vol. 33, No. 3, Dec. 
1979, pp. 405-410. 


8 




Figure 1 . — Rotor 67 test hardware. Figure 3. — Grid leading edge enlargement. 



Figure 2. — Computational grid. 185 C x 40 tangential x 49 
spanwise points. 
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Figure 4. — Comparison of computed and meas- 
ured adiabatic efficiency and total pressure 
ratio characteristic at 100 percent speed. 
Points labeled "near stall" and "near peak 
efficiency" are discussed in detail in the text. 
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Figure 5. — Contours of relative Mach number, near peak efficiency. In computed plots heavy line 
is M = 1 .0, and contour increment is 0.05. 
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Figure 6. — Comparison of computed and measured blade-to-blade distributions of relative Mach 
number, near peak efficiency. 
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Figure 7. — Comparison of computed and measured spanwise distributions of total pressure, total temperature, 
static pressure, and flow angle downstream of the rotor, near peak efficiency. 
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Figure 8. — Computed particle traces, near peak efficiency. 
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Figure 9. — Contours of relative Mach number, near stall. 
1 .0, and contour increment is 0.05. 
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Figure 10.— Comparison of computed and measured blade-to-blade distributions of relative Mach 
number, near stall. 
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Figure 11. — Comparison of computed and measured spanwise distributions of total pressure, total tempera- 
ture, static pressure, and flow angle downstream of the rotor, near stall. 
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(c) Trailing edge 
enlargement, looking 
upstream. 


Figure 12. — Computed particle traces, near stall. 
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